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ABSTRACT 

We investigate the formation of a galaxy reaching a virial mass of « 10® Mq at 2 « 10 
by carrying out a zoomed radiation-hydrodynamical cosmological simulation. This simula¬ 
tion traces Population III (Pop III) star formation, characterized by a modestly top-heavy ini¬ 
tial mass function (IMF), and considers stellar feedback such as photoionization heating from 
Pop III and Population II (Pop II) stars, mechanical and chemical feedback from supernovae 
(SNe), and X-ray feedback from accreting black holes (BHs) and high-mass X-ray binaries 
(HMXBs). We self-consistently impose a transition in star formation mode from top-heavy 
Pop III to low-mass Pop II, and find that the star formation rate in the computational box is 
dominated by Pop III until z ~ 13, and by Pop II thereafter. The simulated galaxy experi¬ 
ences bursty star formation, with a substantially reduced gas content due to photoionization 
heating from Pop III and Pop II stars, together with SN feedback. All the gas within the sim¬ 
ulated galaxy is metal-enriched above 10“® Zq, such that there are no remaining pockets of 
primordial gas. The simulated galaxy has an estimated observed flux of ~ 10“® nJy, which 
is too low to be detected by the James Webb Space Telescope (JWST) without strong lensing 
amplification. We also show that our simulated galaxy is similar in terms of stellar mass to 
Segue 2, the least luminous dwarf known in the Local Group. 

Key words: cosmology: theory - galaxies: formation - galaxies: high-redshift - HII regions 
- hydrodynamics - intergalactic medium - supemovae: physics. 


1 INTRODUCTION 

An important goal of modem cosmology is to understand the for¬ 
mation and properties of early galaxies that formed a few hun¬ 
dred million years after the Big Bang (Bromm & Yoshida 2011; 
Loeb & Furlanetto 2013; Wiklind et al. 2013). In the context 
of hierarchical cosmology, where small stmctures form first and 
evolve into bigger systems, dwarf galaxies with virial masses of 
Mvir = 10^ —1O®M0 ata ~ 6—15 were the basic building blocks 
of large galaxies seen today. Furthermore, it is thought that these 
galaxies played a pivotal role in reionizing the Universe, the global 
phase transition from an early neutral medium to a fully ionized 
one (for reviews see, e.g. Barkana & Loeb 2007; Furlanetto et al. 
2006; Robertson et al. 2010). Specifically, the first galaxies were 
the drivers of the initial stages of reionization (Haardt & Madau 
2012; Shull et al. 2012; Robertson et al. 2013). 

Over the past decade, there have been enormous efforts to 
push the high redshift frontier to 2 ~ 7, detecting light from galax¬ 
ies with a stellar mass of M* > 10® M© (e.g. Bouwens et al. 2011; 
Oesch et al. 2012; Finkelstein 2013; Schenker et al. 2014; Stark 
et al. 2015). To assess the observability of even smaller galaxies, 
and galaxies at even higher redshifts, we need to predict the prop- 
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erties of the first galaxies by pushing numerical simulations to new 
levels of physical realism and detail. Our goal here is to derive ab- 
initio predictions for key quantities, such as their stellar popula¬ 
tion mix, quantified by a possibly time-variable initial mass func¬ 
tion (IMF), star formation rates (SFRs), metallicities, and result¬ 
ing broad-band colour and recombination-line spectra (e.g. John¬ 
son et al. 2009; Pawlik et al. 2011; Zackrisson et al. 2011). A SFR 
above 0.1 M© yr“^, corresponding to M* ~ 10^ M©, might be 
required for galaxies at 2 > 10 to be observed with the James Webb 
Space Telescope (JWST) (e.g. Inoue 2011; Pawlik et al. 2013). For 
JWST deep fields to reach even fainter systems, a magnification 
boost from gravitational lensing may be required (Zackrisson et al. 
2012, 2015). 

The properties of such early galaxies were determined by the 
feedback from preceding generations of stars (for a review see, 
e.g. Ciardi & Ferrara 2005). The first generation of stars, the so- 
called Population III (Pop III), formed at 2 < 30 in dark mat¬ 
ter (DM) minihaloes with Mvir = 10® — 10® M©, predominately 
via ro-vibrationally excited molecular hydrogen (H 2 ) cooling (e.g. 
Haiman et al. 1996; Tegmark et al. 1997; Bromm et al. 2002; 
Yoshida et al. 2003). Radiation from Pop III stars dramatically al¬ 
tered the gas within their host minihaloes, through photoionization, 
photoheating, and photoevaporation (e.g. Kitayama et al. 2004; 
Whalen et al. 2004; Alvarez et al. 2006; Yoshida et al. 2007; Wise 
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& Abel 2008). Once a Pop III star died as a supernova (SN), heavy 
elements were dispersed, enriching the intra-host and intergalac- 
tic medium (IGM), thus initiating the prolonged process of chemi¬ 
cal evolution (reviewed in Karlsson et al. 2013). The resulting SN 
feedback from the first stars greatly affected the gas in the shallow 
potential wells of the host systems (e.g. Greif et al. 2007, 2010; 
Whalen et al. 2008). The more distant, diffuse IGM was heated as 
well by X-rays emitted by accreting single black holes (BHs), or 
high-mass X-ray binaries (HMXBs), both remnants of Pop III stars 
(e.g. Glover & Brand 2003; Kuhlen & Madau 2005; Milosavljevic 
et al. 2009a,b; Alvarez et al. 2009; Mirabel et al. 2011; Wheeler 
& Johnson 2011; Hummel et al. 2014; Jeon et al. 2012, 2014a; Xu 
et al. 2014). 

The critical role of feedback from stars in massive galaxies 
with Mvir ~ 10® — 10^® Mq is also emphasized by many authors 
(e.g. Schaye et al. 2010; Stinson et al. 2013; Kim et al. 2013; Wise 
et al. 2014). Photoionization, radiation pressure, and SN feedback 
operate in a non-linear, coupled fashion, such that considering a full 
treatment of the radiation-hydrodynamics and time dependence of 
these feedback mechanisms is necessary to reproduce the observed 
histories of star formation, including their stochastic nature, par¬ 
ticularly in dwarf galaxies (e.g. Pawlik et al. 2009; Hopkins et al. 
2014). However, due to the prohibitive computational cost related 
to the baryonic physics, simulating the assembly of massive galax¬ 
ies from first principles is not feasible. The first dwarf galaxies with 
virial masses Mvir < 10®, where such ab-initio modeling is com¬ 
ing within reach, thus provide us with ideal laboratories to develop 
a detailed theory of feedback-regulated galaxy formation. 

The birth of the first galaxies is often defined by theorists 
as the emergence of atomic cooling haloes, systems assembled at 
z < 15 with masses of Mvir ^ 5 x Mq, within which the 
primordial gas was able to cool via atomic hydrogen lines instead 
of molecular hydrogen (e.g. Oh c& Haiman 2002; Bromm & Loeb 
2003a). The deeper potential wells of these haloes allowed them 
to sustain self-regulated star formation in a multi-phase medium by 
holding on to the gas driven out by stellar feedback. Another feature 
that renders the atomic cooling haloes viable first galaxy candidates 
is the presence of supersonic turbulence, which is believed to play 
an important role in present-day star formation (McKee & Ostriker 
2007). Such turbulence is expected to develop during the build-up 
of the first galaxies due to the inflow of cold gas along cosmic fila¬ 
ments (e.g. Wise & Abel 2007; Greif et al. 2008; Safranek-Shrader 
et al. 2012). 

The relative importance of Pop III stellar feedback on first 
galaxy formation is determined by the mass range of Pop III 
stars (e.g. Bromm 2013). Stars with initial masses in the range 
of 40 Mq — 140 Mq, or larger than 260 M©, are expected to 
directly collapse into BHs, whereas they undergo pair-instability 
supernovae (PISNe) in the range of 140 M© — 260 M© (Heger 
& Woosley 2002). The PISN range can be extended down to 
~ 85 M© if the progenitor is rapidly rotating (Chatzopoulos & 
Wheeler 2012; Yoon et al. 2012). Finally, stars with masses be¬ 
tween 8 Mq and 40 M© likely died as conventional core-collapse 
SNe (CCSNe), or energetic hypernovae in case of high rotation 
rates (e.g. Umeda & Nomoto 2005). 

Earlier studies had predicted that the first stars formed in iso¬ 
lation, with a characteristic mass of ~ 100 M© (e.g. Abel et al. 
2002; Bromm et al. 2002; Yoshida et al. 2006). However, recent 
high-resolution simulations, exploring the high densities encoun¬ 
tered in protostellar disks, have suggested that the primordial cloud 
can further fragment, giving rise to the possibility of binary or mul¬ 
tiple systems (Turk et al. 2009; Stacy et al. 2010; Clark et al. 2011; 


Prieto et al. 2011; Smith et al. 2011; Greif et al. 2011, 2012; Stacy 
et al. 2012; Dopcke et al. 2013). Radiation-hydrodynamical simu¬ 
lations that take into account the feedback from the protostellar ra¬ 
diation suggest that the first stars typically reached masses of a few 
~ 10 M© (e.g. Hosokawa et al. 2011; Stacy et al. 2012; Hirano 
et al. 2014; Susa et al. 2014). If this is indeed the case, the feed¬ 
back from primordial stars would have been less disruptive than 
previously thought (Kitayama & Yoshida 2005; Whalen et al. 2008; 
Chiaki et al. 2013). Specifically, PISN events would have been less 
frequent, with the majority of Pop III stars now expected to die as 
CCSNe. It is therefore important to revisit the assembly process of 
the first galaxies, taking this change in characteristic Pop III mass 
scale into account. 

Most previous studies have been focusing on the feedback 
from massive Pop III stars (> 100 Mq) on first galaxy formation 
(e.g. Wise & Abel 2008; Greif et al. 2010; Wise et al. 2012; Paw¬ 
lik et al. 2013; Muratov et al. 2013; Xu et al. 2013; Wise et al. 
2014). One prominent result is that even a single PISN with an ex¬ 
plosion energy of FisN ~ 10®® erg, and a metal yield of y ~ 0.5, 
has a dramatic impact on subsequent evolution, evacuating the gas 
within the host halo, and preventing further star formation for about 
a Hubble time, ~300Myr at 2 ~ 10. As the ejected metals fall 
back into the emerging galaxy, the central gas is enriched up to 
Z > 10"® ® Z© (e.g. Wise & Abel 2008; Greif et al. 2010). On 
the other hand, the recovery from the feedback of less massive 
Pop III stars (< 40 M©), followed by CCSNe with Ssn ~ 10®^ 
and y = 0.05, was rather prompt, with recovery timescales of only 
a few ~ 10 Myr (Jeon et al. 2014b). Such prompt recovery allowed 
the first galaxies to experience multiple episodes of star formation 
prior to their assembly (Ritter et al. 2012, 2014). Related to this 
weakened feedback, the star-forming gas inside the first galaxies 
was possibly enriched by multiple SNe. The resulting metal enrich¬ 
ment was sufficient to affect a transition of the star formation mode 
from Pop III to normal. Population II (Pop II) stars. Thus, regions 
that hosted moderate mass Pop III stars would locally undergo the 
transition to (low-mass dominated) Pop II star formation earlier on. 

The minimum (“critical”) metallicity for the Pop III/Pop II 
transition depends on which cooling is responsible for the gas frag¬ 
mentation (e.g. Omukai 2000; Bromm et al. 2001; Schneider et al. 
2002; Bromm & Loeb 2003b; Omukai et al. 2005; Schneider & 
Omukai 2010; Schneider et al. 2012; Chiaki et al. 2014; Ji et al. 
2014). Fine-structure line cooling, mainly via C II and O I, above 
Z ~ 10“®'® Z©, leads to vigorous fragmentation of the gas cloud, 
giving rise to the formation of stellar clusters. Further fragmen¬ 
tation can be achieved at very high densities, nn ^ 10^® cm“®, 
via dust continuum cooling if the initial gas-to-dust ratio exceeds 
JJcrit = [2.6 — 6.3] X 10“® (Schneider et al. 2012; Chiaki et al. 
2014). Globally, the Pop III/Pop II transition is a gradual process, 
as the ejected metals are brought back into virialized host haloes, 
keeping a large fraction of the IGM still pristine. Therefore, those 
two populations co-existed, such that Pop III star formation was 
still possible even at 2 ~ 6 (e.g. Trent! et al. 2009; Maio et al. 
2010, 2011; Scannapieco et al. 2003; Wise et al. 2012; Johnson 
etal. 2013). 

Here, we present the results of a cosmological zoom-in, high- 
resolution radiation-hydrodynamics simulation, following the first 
galaxy assembly process from first principles. We consider real¬ 
istic descriptions of Pop III/Pop II star formation, photoionization 
and photoheating from stars, the mechanical and chemical feedback 
from PISNe and CCSNe, and the X-ray feedback from accreting 
BHs and HMXBs. Guided by the results from recent investigations 
of primordial star formation, we in particular explore the implica- 
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tions of a somewhat less top-heavy Pop III IMF, with a character¬ 
istic mass of a few ~ 10 Mq, giving rise to a large fraction of 
CCSNe rather than PISNe. 

The outline of the paper is as follows. Our numerical method¬ 
ology is described in Section 2, and the simulation results are pre¬ 
sented in Section 3. Finally, our main findings are summarized in 
Section 4. For consistency, all distances are expressed in physical 
(proper) units unless noted otherwise. 


2 NUMERICAL METHODOLOGY 

2.1 Gravity, hydrodynamics, and chemistry 

We have performed our simulations using a modified version of the 
W-body/TreePM Smoothed Particle Hydrodynamics (SPH) code 
GADGET (Springel et al. 2001; Springel 2005). The initial con¬ 
ditions within a cubic volume 1 comoving Mpc on a side are gen¬ 
erated by assuming a ACDM cosmology with a matter density pa¬ 
rameter of Dm = 1 — Ga = 0.3, baryon density Gb = 0.04, 
present-day Hubble expansion rate Hq = 70 km s“^Mpc“^, spec¬ 
tral index Ua = 1.0, and normalization as = 0.9. The initial con¬ 
ditions are hierarchically refined by a consecutive zoom-in tech¬ 
nique, rendering the masses of dark matter (DM) and SPH particles 
in the highest resolution region with an approximate linear size of 
300 comoving kpc moM ~ 33 Mq and mspH ~ SMq, respec¬ 
tively. We adopt the gravitational softening length Cgoft = 70 co¬ 
moving pc for both DM and baryonic particles. We start our simu¬ 
lations at z ~ 100 and follow the cosmic evolution of the DM and 
gas until the assembly process of the first galaxy has been com¬ 
pleted at z « 10. 

We consider all relevant primordial chemistry and cooling 
processes. Specifically, we include H and He collisional ioniza¬ 
tion, excitation and recombination cooling, bremsstrahlung, inverse 
Compton cooling, and collisional excitation cooling of H 2 and 
HD. The code self-consistently solves the non-equilibrium chem¬ 
istry for 9 abundances (H, H~, H 2 , HJ, He, He"*", and 

e“), as well as for the three deuterium species D, D+, and HD. 
Additionally, metal-enriched gas can also be cooled by metal 
species, specifically C, O, and Si, with solar relative abundances. 
For the metal cooling rates, we adopt the results of Glover & 
Jappsen (2007). The chemical network comprises the key species, 
C, C^, O, O^, Si, Si"*" and Si^”*". We assume that all the metals 
are in gas phase atoms or ions, and thus we do not consider dust 
cooling as the gas densities explored in this work are far below the 
threshold value above which it starts to become important. 


2.2 Star formation physics 

2.2.1 Population III 


For the star formation, we do not form stars by following the evo¬ 
lution from the initial protostars to the final mass (e.g. Hirano et al. 
2014), rather we adopt a sink algorithm (Johnson & Bromm 2006), 
where stars are represented by sink particles. Once a gas particle 
exceeds the threshold density, nn.max = 10'* cm“®, the highest- 
density SPH particle is converted into a collisionless sink particle, 
subsequently accreting neighboring gas particles until the mass of 
the Pop III star is reached. The masses of Pop III stars are randomly 
sampled from an IMF with a functional form of 


dN 

dlnM 


oc M exp 


( McharV ' 

V M J 


( 1 ) 


where Mchar = 20 Mq is the characteristic mass. Above Mchar, 
it behaves as a Salpeter-like IMF, but is exponentially cutoff below 
that mass (e.g. Chabrier 2003; Wise et al. 2012). Once a sink par¬ 
ticle forms, we prevent the subsequent star formation in a sphere 
with radius 30 pc centered on the sink particle for ~ 1 Myr. This is 
because it takes time for the hydrodynamical shocks with a velocity 
~30 km/s, generated by photoheating to sweep up and evacuate the 
central gas, reducing the central gas densities. Without such tempo¬ 
rary suppression in star formation, multiple stars might form at the 
same location at the same time. This would correspond to frag¬ 
mentation within the SPH smoothing kernel, which is not resolved, 
and would thus be physically unreliably. We allow sink particles to 
merge once their distance falls below 1 pc, similar to the baryonic 
resolution scale (res [(d-YALes}/{47mri, max « 0.5 pc, 

where X = 0.76 is the hydrogen mass fraction. Here, the bary¬ 
onic mass resolution is M^es = (VngbmsPH ~ 240 Mq, where 
Angb = 48 is the number of particles in the SPH smoothing ker¬ 
nel (Bate & Burkert 1997). The new positions and velocities of the 
merged sinks are computed by mass-weighted averaging. 


2.2.2 Population II 

Pop II stars are formed out of the gas cloud which is metal-enriched 
by the previous generation of stars. The star formation recipe is 
similar to that of Pop III star formation except for one more con¬ 
dition: if the metallicity of a gas particle, eligible for star forma¬ 
tion, exceeds the characteristic metallicity, Zcv\t = 10“^'® Zq, 
for the transition from Pop III to Pop II, we assume that a newly 
formed sink particle represents a Pop II star cluster. Therefore, once 
a SPH particle satisfies those two criteria, nH,max = 10* cm“® and 
Zcvit = 10“^'® Zq, we immediately create an effective sink parti¬ 
cle with a mass of Mpop ii = 500 M©, by accreting surrounding 
gas particles. The mass of a Pop II cluster Mpop n = 500 M© 
is chosen to host one massive star ~ 20 M©, given a Salpeter 
IMF, dN/dlogm « m~°‘, with a slope a = 1.35, which is from 
Schaerer (2003) where they provide photoionization rates derived 
from an assumed IMF for a Pop II cluster. In the presence of a 
20 Mq star, further star formation inside of the cluster after the 
initial burst is truncated by the photoheating of the star (Safranek- 
Shrader et al. 2014b). 


2.3 Feedback physics 

2.3.1 Photoionization feedback 

Once a star forms, the star particle emits ionizing radiation over 
its lifetime with a black-body spectrum. We use the radiative trans¬ 
fer (RT) code TRAPHIC to transport ionizing photons (Pawlik & 
Schaye 2008, 2011). In simulations with TRAPHIC, photon pack¬ 
ets from radiation sources are transferred along the irregular, spa¬ 
tially adaptive grid set by SPH particles in a photon-conserving 
manner through the simulation box. Thanks to a photon packet 
merging technique in TRAPHIC, the computational cost is inde¬ 
pendent of the number of ionizing radiation sources, which enables 
large simulations considering many sources, and as presented here. 
The RT is controlled by a set of parameters, which we choose iden¬ 
tical to those employed in Jeon et al. (2014a), except that we here 
use two frequency bins, = 2 instead of = 4, with bounding 
energies located at [13.6 eV, 400 eV, 10 keV]. While the main fo¬ 
cus of Jeon et al. (2014a) is to investigate the impact of X-rays on 
the gas in haloes and IGM, requiring to track high energy photons 
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denoted by the high characteristic frequencies, in this work we fo¬ 
cus on the stellar feedback impact. To do so, two frequency bins, 
Ni, = 2, are sufficient to represent the stellar spectrum. We refer 
the reader to the original papers (Pawlik & Schaye 2008, 2011) for 
further details concerning RT methods and to Jeon et al. (2014a) 
for the RT parameters we use in this work. The RT computation is 
coupled to the hydrodynamical evolution by passing the photoion¬ 
ization, photoheating, and photodissociation rates computed by the 
RT to the non-equilibrium solver for the chemical and thermal evo¬ 
lution of the gas. 

For the ionizing photon rates of Pop III stars and their life¬ 
times, we employ polynomial fits (Schaerer 2003) as a function 
of initial mass of a star, log^Q Nio-n ~ 43.61 -I- 4.90* — 0.83*^ 
and logjof* = 9.785 — 3.759* -I- 1.413*^ — 0.186*^, respec¬ 
tively, where * = logjQ(Mpop iii/ Mq). For Pop II clusters, we 
use the value from Schaerer (2003) that provides the properties 
of integrated stellar populations at various metallicities. In partic¬ 
ular, we fixed the ionizing photon rate for the initial ~ 3 Myr 
at (Vion = 10 '^^s“^ M 0 ~^, which is derived by averaging over 
a Salpeter IMF with a slope of a = 1.35 and a metallicity of 
Z = 5 X 10“® Zq in the mass range of [1 Mq, 150 M©]. For the 
following evolution of the rates, we use a polynomial fit (Schaerer 
2003), logjg iVion = 48.13 — 0.42f + O.Oli^ M©”^, where t 
is the elapsed time in Myr after Pop II cluster formation. The emis¬ 
sion from Pop II stars lasts for 20 Myr, corresponding to the typical 
lifetime of OB stars in the cluster. 

In addition to the ionizing radiation from individual Pop III 
stars and Pop II clusters, we also track the Lyman-Werner (LW) ra¬ 
diation from both populations that dissociates molecular hydrogen 
(H 2 ) and deuterated hydrogen (HD). This is treated in the opti¬ 
cal thin limit with a self-shielding correction (Wolcott-Green et al. 
2011). The rate at which LW photons are emitted by Pop III stars is 
fit using logjg JVlw = 44.03 -I- 4.59* — 0.77*^ (Schaerer 2002). 
We only allow stars to form within the highly resolved region corre¬ 
sponding to a linear size of « 300 comoving kpc and do not follow 
the propagation of photons outside the region. For simplicity, we 
do not consider a global UV background here. 

2.3.2 Supernova feedback 

The SN explosion energy is inserted as thermal energy into the 
neighboring SPH particles, iVneigh = 32, a subset of the parti¬ 
cles within the SPH smoothing kernel, A'^neigh, around a sink par¬ 
ticle. To ensure high accuracy and energy conservation, we use a 
timestep-limiter, with which the ratio of timesteps of neighboring 
SPH particles cannot be larger than 4 (Saitoh & Makino 2009, 
Durier & Dalla Vecchia 2012). Additionally, using a timestep- 
updater, implemented hy Durier & Dalla Vecchia (2012), all neigh¬ 
boring particles, iVneigh, around an exploding sink particle become 
active particles at the time of SN energy injection, enabling an im¬ 
mediate hydrodynamical response. 

At the end of a Pop III star’s lifetime, the fate of the Pop III 
star is determined by its initial mass. For instance. Pop III stars with 
masses between 9 M© and 40 M© are expected to end their lives 
as conventional core collapse SNe. Here, we use the fixed CCSN 
energy of iSsN = 10®^ ergs and metal yield, y = 0.05, for Pop III 
progenitors. We should mention that the CCSN explosion energies 
and the metal yields can vary depending on the progenitor’s mass 
(e.g. Heger & Woosley 2010). In particular, under the existence of 
rotation of Pop III progenitor in a mass range of 25 — 40 M©, a 
much stronger explosion may be expected, likely resulting in a hy¬ 
pernova (e.g. Fryer & Heger 2000; Nomoto et al. 2006). However, 
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Figure 1. The history of comoving star formation rate density (SFRD) in 
the simulated region for Pop III (black solid line) and Pop II (red dotted line) 
stars. For comparison, we also show the analytic fitting formula, derived by 
Hernquist & Springel (2003), for higher-mass systems where atomic cool¬ 
ing is active (light blue-dashed line). Earlier on, at z > 13, Pop III star 
formation dominates, but as the gas is metal-enriched, the Pop II mode be¬ 
comes dominant below 2 : r; 13. 


due to the uncertainty in the degree of spin of Pop III stars, which is 
subject to stellar winds, gravitational and hydromagnetic instabili¬ 
ties (e.g. Stacy et al. 2011), we only consider conventional CCSNe 
by fixing the explosion energy and metal yield. 

For the highly energetic PISNe from Pop III stars with masses 
in the range of 140 M© < M* < 260 M© (Heger & Woosley 
2002), we adopt the SN energy of 10®^ ergs and y = 0.5. Out of 
a Pop II cluster with a mass of 500 Mq, a total of 9 stars more 
massive than 8 Mq , above which the stars are eligible to explode 
as CCSNe, are expected to form. We simultaneously trigger 9 SNe 
explosions at the end of the lifetime of a Pop II cluster, injecting an 
energy of 9 x 10®^ ergs with an IMF-averaged effective metal yield 
yeff = 0.005 onto the neighboring particles around the Pop II sink 
particle. 

Initially, the metals from SN explosions are evenly distributed 
among the neighboring SPH particles, implying initial metallicities, 

rTT-metal.i 

msPH + rUmetal,! 

where mmetau = M^y/Nneigh for Pop III stars and mmetai.i = 
Mpop iiyeff /Nneigh for Pop II clusters. The lack of mass flux be¬ 
tween SPH particles means that metals are transported ballistically 
in an SPH scheme (e.g. Wiersma et al. 2009). Physically, the mix¬ 
ing of metals is achieved through an unresolved turbulent cascade 
from large scales to small scales, which is unresolved. Therefore, 
often, a metal mixing is implemented by subgrid models in which 
physical processes below the smallest resolved scale are deter¬ 
mined by physical quantities at resolved scales (e.g. Schmidt & 
Federrath 2011). We adopt a diffusion-based method implemented 
by Greif et al. (2009) which assume that the mixing efficiency on 
unresolved scales is determined by the physical properties on the 
scale of the SPH smoothing kernel (Klessen & Lin 2003). Conse¬ 
quently, the metal diffusion rate is set by the local diffusion coeffi¬ 
cient, D = 2 p V I, where the length scale, I, is comparable to the 
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smoothing length of the SPH kernel, I — h, and p is the gas density. 
The velocity dispersion within the kernel, v, is given by 


-2 


Vi 


1 

A^ngb 


'^\Vi-Vj\^. 


(3) 


Here, Vi and Vj are the velocities of particles i and j within the 
kernel. 


2.3.3 Black hole feedback 


Pop III stars with masses between AO Mq and 140 M©, or larger 
than 260 M© will directly collapse into BHs (Heger et al. 2003). 
The BHs then grow by Bondi-Hoyle accretion (Bondi & Hoyle 
1944), assuming that a fraction e = 0.1 of the accreted mass is 
converted into ionizing radiation. Such stellar mass BHs accreting 
diffuse gas are called miniquasars (Kuhlen & Madau 2005). The 
ionizing luminosities of miniquasars are determined by normaliz¬ 
ing the total accretion luminosity 


^10keV//ip 

Tbh = / Ljjdv = - -Mbhc^, (4) 

Jo 1 - e 

where hp is Planck’s constant and c is the speed of light. The BH 
(sink) particles grow in mass by accreting surrounding gas accord¬ 
ing to the accretion rate, Mbh = (1 — e)Macc, which is computed 
as 


TttG^MIhP 

(ci+u2J3/2- 


(5) 


Here, Cs is the sound speed, p the gas density, Mbh the BH mass, 
and Wrei the relative velocity of the BH to the surrounding gas. Note 
that Bondi-Hoyle accretion is thus a numerical input, or recipe, and 
not an independent simulation result. The estimated rate is an up¬ 
per limit. In reality, the true accretion rates are likely smaller than 
the nominal Bondi-Hoyle value, if feedback from radiation pres¬ 
sure is considered (e.g. Milosavljevic et al. 2009a,b; Park & Ricotti 
2012). In addition, at later time as the gas in an atomic cooling 
halo becomes turbulent, the vorticity of the turbulent gas can re¬ 
duce the accretion onto the BH (Krumholz et al. 2005). Motivated 
by observations at low redshifts, we take the same form of a spec¬ 
tral energy distribution of the emerging radiation from an accreting 
single black hole, which is characterized by a thermal multi-color 
disk (e.g. Pringle 1981; Mitsuda et al. 1984) at frequencies lower 
than 0.2 keV /hp and a non-thermal power law component at higher 
frequencies (e.g. Kuhlen & Madau 2005). For the details, we refer 
the reader to Jeon et al. (2014a). 

Further, we assume that half of the black holes are formed in 
a binary system and then every third of those evolve into a HMXB 
after the primary turned into a BH (Power et al. 2009). For the 
HMXB phase of duration AIhmxb = 2 Myr, corresponding to the 
typical main-sequence lifetime of the donor stellar companion (e.g. 
Belczynski et al. 2012), we assume an HMXB luminosity equal to 
the Eddington luminosity. 


nlOkeV/tip 


Thmxb = 


1.4 X 10‘‘°erg 


Lirdiy — T/Edd 

Mbh 


100 M© J ’ 


( 6 ) 

(7) 


corresponding to accretion of gas from the companion at a rate 
2.2 X 10“® M© yr~^(MBH/100M©). The X-ray radiation is fol¬ 
lowed using RT. We also take into account the secondary ioniza¬ 
tions, which are mainly done by high energy electrons left behind 


the first ionizations of hydrogens by hard energy photons emitted 
from an accreting BH or a HMXB. We use the energy dependent 
fits to the secondary ionization and heating fractions provided by 
Ricotti et al. (2002). 


3 SIMULATION RESULTS 

In the following, we present our results. First, in Section 3.1, we 
discuss the global properties of the IGM and of haloes within the 
highly resolved region and their evolution to 2 ; « 10. In Section 
3.2, we focus on the assembly process of the central halo hosting 
the emerging galaxy. Then, based on the physical quantities derived 
from this galaxy, we compare the properties of the first galaxy with 
those of the local counterparts in Section 3.3. Finally, in Section 
3.4, we discuss the detectability of the simulated first galaxy. 


3.1 Global evolution 

3.1.1 Star formation history 

In the simulation, the first Pop III star, with a randomly assigned 
mass of 100 M©, forms at 2 « 28 out of the primordial gas in a 
minihalo with mass Mvir ~ 5 x 10® M©. After its short lifetime 
of 2.7 Myr, the star ends up as a BH in a binary system, forming a 
HMXB. Over a brief interval of AIhmxb = 2 Myr, X-rays emit¬ 
ted from the HMXB further heat the surrounding gas, which was 
already heated to ~ lO"* K by the Pop III progenitor, and they heat 
the distant IGM to a few ~ 10® K. At 2 « 23, the first Pop II star 
formation event takes place in metal-enriched gas, which was pre¬ 
viously polluted by a Pop III star, whose mass of 20 M© was again 
randomly drawn from the underlying IMF Pop III star. The time de¬ 
lay between the onset of second-generation. Pop II, star formation 
and the death of the 20 M© Pop III star is ~ 16 Myr, consistent 
with our earlier study (Jeon et al. 2014b). By the end of the sim¬ 
ulation at 2 « 10.5, a total of Npop n = 65 Pop II clusters and 
Npop III = 183 individual Pop III stars has formed, with 128, 2, 
44, and 11 Pop III stars ending their lives as CCSNe, PISNe, BHs, 
and HMXBs, respectively. 

Fig. 1 shows the comoving star formation rate density (SFRD) 
for the Pop III individual stars and the Pop II clusters. During the 
initial 50 Myr, only Pop III stars exist. However, once metal pollu¬ 
tion commences by Pop III SNe, the metallicity of the gas increases 
to trigger the transition of the star formation mode from Pop III to 
Pop II. Therefore, at 2 < 23, two star formation modes coexist over 
the next ~ 350 Myr. As the dense gas available for star formation 
gets metal-enriched, the dominant star formation mode changes 
from Pop III to Pop II. Accordingly, the SFRD for Pop II stars 
begins to increase and reaches a few ~ 10“® M©yr“®Mpc“® 
at 2 « 14. On the contrary. Pop III star formation starts to de¬ 
cline from ~ 7 X 10“® M©yr“®Mpc“® at 2 « 14 to a few 
~ 10“® M©yr“^Mpc“® at 2 « 10.5 as the amount of dense pris¬ 
tine gas decreases. 

It is noteworthy that the SFRDs measured in the simulation 
are very episodic rather than continuous. This episodic, bursty star 
formation is mainly caused by a stellar feedback cycle: dense gas 
inside a halo is evacuated by photoheating from stars, the subse¬ 
quent SN blastwave, or both, rendering it gas-starved for a certain 
time, until fresh gas is provided through smooth accretion or merg¬ 
ers with other haloes. Such self-regulated star formation on short 
timescales is also exhibited by simulations of more massive galaxy 
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Figure 2. Morphology in the highly resolved region at 2 Ri 10.5. Clockwise from upper left. Dark matter overdensity, H II fraction, gas metallicity, and gas 
temperature, averaged along the line of sight within the central ~ 350 kpc (comoving). The central white circles denote the size of the virial radius, Ri 1.1 
kpc, of the emerging target galaxy. The target halo grows in mass through mergers with neighboring haloes and accretion along the filaments of the cosmic 
web. The gas in haloes that are undergoing star formation is evacuated and heated, and metal-enriched by stellar feedback. The main metal contributors to the 
diffuse IGM are two energetic PISNe, formed at 2 R 16 and 2 R 15.5, and ejecting a total of ~ 160 Mq in heavy elements into the IGM. 


formation (Hopkins et al. 2014), where stellar feedback, i.e. photo¬ 
heating, SNe, winds, and radiation pressure, is taken into account 
with approximate sub-grid prescriptions. For comparison, we over¬ 
plot the analytic fitting formula, derived by Hernquist & Springel 
(2003), for higher-mass systems, where atomic cooling is active, 
which agrees well with the SFRD for Pop II stars predicted by 
our simulation. We, however, should note that the SFRDs both for 
Pop III and Pop II stars are higher than those computed in the adap¬ 
tive mesh refinement simulation of Wise et al. (2012) by a factor 
of ~ 5. The reason for this is the biased nature of our zoomed, 
high-resolution region, focusing on a small volume and omitting 
under-dense void regions, which results in enhanced star formation 
activity compared to the cosmic mean. 


3.1.2 Global impact of stellar feedback 

Fig. 2 shows images of the dark matter overdensity, ionized hydro¬ 
gen fraction, gas temperature, and gas metallicity at 2 R 10.5 in 
the high-resolution region, corresponding to ~ 350 comoving kpc 
across, centered on the site of the emerging target galaxy. The virial 
radius of the halo hosting the target galaxy, depicted as the white 
dashed circle in Fig. 2, grows over time from ~ 100 pc at 


2 = 28 to Tvir ~ 1.1 kpc at 2 = 10.5. Correspondingly, the halo 
grows in mass from 5 x 10® Mq to 8 x 10^ Mq. Both the gas in¬ 
side haloes and in the diffuse IGM are subject to a variety of heating 
sources, such as UV radiation from Pop III and Pop II stars. X-rays 
from accreting BHs and HMXBs, and thermalized SN explosion 
energy. The temperature inside HII regions reaches T r 10'* K by 
stellar photoionization heating, and further increases to 10^ — 10* 
K due to SN shock heating. Such high-temperature gas is primar¬ 
ily cooled via inverse Compton cooling off the cosmic microwave 
background (CMB) and free-free emission to 10® K, and then H 
and He recombination cooling begins to dominate at the lower tem¬ 
perature of 10"^ - 10® K (Greif et al. 2007). 

In order to disentangle the contribution of radiation sources to 
ionization and heating, in Fig. 3 we compare the ionizing photon 
rates from Pop III stars. Pop II stars, and HMXBs. While the ioniz¬ 
ing photons from Pop III stars are responsible for ionizing the gas at 
high redshifts 2 > 13, the contribution from Pop II stars becomes 
competitive at later times. On the other hand, the contribution from 
a handful of HMXBs is negligible owing to their rarity and their 
short lifetime, AIhmxb = 2 Myr, although the latter quantity is 
very uncertain and may be larger. In Jeon et al. (2014a), the impact 
of X-rays from HMXBs on star formation was similarly negligible. 
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Figure 3. The evolution of ionizing photon rates (/eji) and total accumulated photons (right) from the Pop III stars. Pop II clusters, and HMXBs. The 
contribution from Pop III stars dominates the total ionizing photon budget up to 2 Ri 10.5, but Pop II clusters are about to overtake it, as the global star 
formation mode changes to Pop II star formation. The total photon production from the 11 HMXBs formed during the run is negligible, mainly due to our 
assumption here of a very brief interval of 2 Myr, over which the X-ray source is active. This interval, however, is uncertain, and could be larger. If so, the 
relative HMXB contribution to the ionizing photon budget could be more significant. 


while the effect of X-ray heating in the diffuse IGM was signifi¬ 
cant, smoothing out small-scale structures. This smoothing effect 
is less prominent in the current simulation, which hosts only a sin¬ 
gle HMXB by 2 ~ 18, instead of the 8 sources present at that time 
in our previous study. 

We assess the impact of the stellar feedback on individual 
haloes by using the SUB FIND halo finder (Springel et al. 2001) to 
identify haloes from the highly-refined region in the simulated vol¬ 
ume, allowing us to investigate their gas, stellar, and metal prop¬ 
erties. Each circle in Fig. 4 marks the mass-weighted average of 
the corresponding stellar or gas properties inside the virial sphere, 
which is centered on the most bound DM particle of the halo, and 
bounded by radii where the DM density is 200 times higher than 
the mean density of the Universe. Star formation takes place only 
in haloes more massive than Mvir = 5 x 10® M©, below which the 
haloes are not massive enough to satisfy the condition fcooi < Iff, 
where tcooi and ta are the cooling time and the free-fall time, re¬ 
spectively, for the gas to cool and form stars (e.g. Tegmark et al. 
1997; Bromm et al. 2002). 

Fig. 4 shows the gas properties of the extracted haloes at 
2 « 10.5, including the baryon fraction, the total mass of Pop III 
BH remnants and Pop II low-mass (long-lived) stars, the highest 
gas density achieved, and the average metallicity of the gas in the 
halo. Photoheating and SNe drive the gas out of the haloes, lower¬ 
ing the baryon fraction below the cosmic mean (upper-left panel of 
Fig. 4). Among the haloes more massive than 10® Mq atz~ 10.5, 
about 40%, corresponding to 7 out of 17 haloes (blue circles in 
Fig. 4), experience photoheating and SN blastwave heating by both 
Pop III stars and Pop II clusters. Specifically, the three most mas¬ 
sive haloes, (Mvir ~ 7 x 10^ Mq, ~ 3 x 10^ Mq, ~ 2 x 10^ M©) 
at 2 « 10.5 have experienced (32, 9, 7) Pop III SNe and (28, 2, 
3) SNe from Pop II clusters. The gas fractions of these haloes at 
2 « 10.5 are 5.1%, 3.3%, and 3%, respectively. For lower-mass 
haloes, on the other hand, the scatter in the gas fraction is much 
more substantial, due to their susceptibility to negative stellar feed¬ 
back in the shallow potential wells. Those haloes are more likely to 


lose their gas once they undergo star formation, and it takes longer 
for the ejected gas to be re-incorporated. 

About 35% of the haloes larger than Mvir ^ 10® M©, de¬ 
noted as red circles, only host Pop III stars. Most of them expe¬ 
rience subsequent CCSN explosions, whereas 2 haloes host a rare 
PISN. The effect of a PISN is much more disruptive than that of 
a CCSN. As an example, the third heaviest halo hosts a 160 M© 
PISN at 2 = 16, when its mass reaches Myir ~ 10^ M©. This 
hyper-energetic explosion effectively evacuates the gas within the 
halo, lowering the baryon fraction by an order of magnitude, from 
/bar = 0.088 at 2 = 16 to /bar = 0.0046 at 2 = 13.7. Subse¬ 
quently, the halo gas begins to be replenished through accretion 
of both the ejected metal-enriched gas and fresh primordial gas 
along the cosmic filaments, such that at the end of the simulation 
the baryon fraction again reaches /bar = 0.03. Despite the gas re¬ 
covery, however, any further star formation inside the halo has been 
completely suppressed since 2 = 16. 

Haloes with Mvir < 10® M©, corresponding to the mass 
range considered in this work, are susceptible not only to the lo¬ 
cal stellar feedback inside the host itself, but also to the external 
global UV background that also photoheats the gas in their shal¬ 
low potential wells, reducing the baryon fraction. For comparison, 
we overplot a fit (dotted line the upper-left panel of Fig. 4) that 
describes how the baryon fraction is affected by reionization as a 
function of halo mass and redshift (Cnedin 2000; Okamoto et al. 
2008): 


/bar(M, 2) =< /bar > 



1 ) 


M 

Me (2) 


- 3/13 


,( 8 ) 


where </bar> is the cosmic mean baryon fraction, /I = 2 a fitting 
parameter, and Me (2 = 10) = 7x 10® M© the characteristic mass, 
below which galaxies are strongly affected by the UV background. 
The gas in haloes affected by internal feedback from both Pop III 
and Pop II stars tends to be more strongly evacuated, compared to 
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Figure 4. The properties of gas, stars, and metals in haloes at 2 ps 10.5, extracted from the high resolution region using the SUBFIND halo finder (Springel 
et al. 2001). Clockwise from upper left: halo baryon fraction, total mass in Pop III BH remnants and long-lived (low-mass) Pop II stars, average metallicity of 
the total halo gas and of the high-density gas, and the maximum hydrogen number density. These properties are computed by mass-weighted averaging inside 
the virialized region centered on the most bound halo particle. Due to the feedback from photoheating, SNe, or both, the gas inside the haloes is significantly 
evacuated, giving rise to baryon fractions below the cosmic mean value (dashed horizontal line in the upper-left panel). For comparison, we also show the 
predicted suppression in baryon fraction due to the presence of a global, external UV background, as given in Equ. 8 (dotted line in the upper-left panel). 
About 46% of star-forming haloes (colored circles) only host a single or a few Pop III stars (red circles), whereas the remainder hosts both populations (blue 
circles). The haloes more massive than Mvir = 10® Mq are substantially metal-enriched, giving rise to an average metallicity above Z = 10“®'® Zq. 
The dotted and dashed horizontal lines in the bottom right panel denote the critical metallicity for the Pop III/Pop II transition, Zcrit = 10-3-5 Z0 and 
■^crit.dust = 10“5'5 Zq, Corresponding to metal fine-structure line cooling and dust continuum cooling, respectively. Depending on the choice of critical 
metallicity, the timing of the Pop III to Pop II transition in a given halo can vary. 


the external UV background case, whereas /bar in haloes hosting 
only Pop III stars is less severely attenuated. 

The average gas metallicity in the haloes with Mvtr ^ 
10® Mq lies between 2avg = 10“® — IO^^Zq. The metal- 
enrichment in the individual haloes is mainly achieved by CCSNe, 
because the relatively low characteristic mass, Mchar = 20 Mq, 
for the Pop III IMF results in only two PISNe during the run. 
When the 160 Mq PISN is triggered, the surrounding gas is im¬ 
mediately polluted to an average of2 = 4xlO~®ZQ. However, 
those heavy elements are quickly driven out beyond the virial ra¬ 
dius of the host halo. After 50 Myr, the ejected metals start to fall 
back as gas is accreted, and the average gas metallicity reaches 
2 = 3 X 10“"^ Zq in the halo with a mass of Mvir ~ 2 x 10^ Mq 
at 2 « 10.5. This agrees with the results by Wise et al. (2012), 
where they found that for haloes Mvir > 3 x 10^ Mq, all haloes 
are enriched above [Z/H] > —4. However, we should point out 
that the metal contributors were Pop III PISNe in their work, while 
the haloes at 2 = 10.5 in the current work are mainly polluted by 
CCSNe. At 2 « 10.5, if we consider all haloes with a mass above 


Mvir = 10® Mq, only 2% of them experience in-situ star forma¬ 
tion, while 25% of the haloes exhibit an average gas metallicity 
above Z = 10“® Zq. Thus, some of them were enriched by SNe 
in nearby haloes. 

Given that star formation takes place in cold dense gas, specif¬ 
ically considering the metallicity of the densest gas in the haloes is 
meaningful. In the most massive three haloes, the densest gas is 
enriched above Zcrit = lO”®'® Zq (see the bottom-right panel of 
Fig. 4), indicating that Pop II stars are expected to form out of the 
gas. However, it does not necessarily imply that the haloes cease 
forming Pop III stars. Instead, Pop III stars may continue to form 
until the end of the simulation in the highest density regions, which 
are already metal-enriched, but possibly insufficiently so to enable 
Pop II star formation. The time at which the Pop III/Pop II transi¬ 
tion occurs depends on the choice of the critical metallicity, which 
is still quite uncertain. For instance, the densest clouds (open circles 
above 10® cm” 3 in the bottom-left panel of Fig. 4) in the haloes 
between 3 x 10® Mq < Mvir < 10^ Mq are about to form stars. 
Based on Zcrit = 10“®'® Zq, which is assumed in this work, some 
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Figure 5. History of the global metal enrichment. As a series of SNe are 
triggered, the total amount of metals in the simulated region increases over 
time. Due to the large yield, y = 0.05, of Pop III CCSNe, compared to the 
effective yield of Pop II SNe, peff = 0.005, the former dominate, a situa¬ 
tion that is reversed at later times, 2 : < 11. In terms of the total metal budget, 
about 75% of the heavy elements are produced by only two PISNe owing 
to their extremely lai'ge yield, y ^ 0.5, in combination with their large pro¬ 
genitor mass. The average metallicity increases up to Zavg = 10“'^ Zq at 
2 : 5 ^ 10.5. On the other hand, the majority of the simulated volume remains 
almost pristine, below Z < 10“^'^ Zq or Z < 10“^'^ Zq, as shown in 
the bottom panel. 


Figure 6 . Top: Evolution of the metal-enriched volume filling fraction as 
a function of metallicity at four different times. As more and more regions 
are polluted by SNe, the volume filling fraction increases over time, but 
the peaks remain at the same metallicity Z = 10“^ — 10“^ Zq at all 
times. We normalize the volume by the total simulated volume, most of 
which is still pristine as shown in the bottom panel of Fig. 5. Bottom: The 
distribution of metallicity depending on the gas density. Which range of gas 
density is preferentially metal-enriched is sensitively determined by the gas 
conditions prior to a SN explosion. For a SN with a low mass progenitor, 
high density gas (tih = 100 cm~^) is likely to be polluted, while the low 
density gas (tih = 0.1 — 0.001 cm~^) is preferentially enriched by a SN 
triggered by a massive progenitor. 


Pop III stars are expected to form, whereas Pop II stars would form 
if we adopted a critical metallicity, Zcrit = 10 -®" Zq (dashed 
horizontal line in the bottom-right panel of Fig. 4), set by dust- 
continuum cooling (e.g. Omukai 2000; Bromm et al. 2001; Schnei¬ 
der et al. 2002; Bromm & Loeb 2003b; Omukai et al. 2005; Schnei¬ 
der & Omukai 2010; Schneider et al. 2012; Chiaki et al. 2014; Ji 
et al. 2014). Therefore, the latter case would give rise to an earlier 
transition in star formation, and would more rapidly extinguish the 
Pop III mode. 


3.1.3 Metallicity evolution 

By the end of the simulation, a total of 128 out of 183 Pop III stars 
end up as CCSNe and 2 stars explode as PISNe. Fig. 5 displays the 
resulting metal enrichment history of the simulated region. While 
the gas confined inside the haloes is mainly enriched by Pop III 
and Pop II CCSNe, the enrichment in the diffuse IGM is domi¬ 
nated by the rare PISN events. Indeed, although only two PISNe 
have occurred, they produce a total of ~ 150 Mq in metals, ac¬ 
counting for 75% of the total metal budget at z ~ 10.5. Intrigu- 
ingly, the optimal strategy to search for PISN-enriched material 
may thus be to probe the high-redshift IGM with absorption spec¬ 
troscopy (Wang et al. 2012). The average mass-weighted metallic¬ 


ity of the region increases to 10“'* Zq at 2 ; ~ 10.5. However, the 
majority of the region is metal-enriched below the critical metal¬ 
licity Zcrit = 10“^ ® Zq, giving rise to a volume filling factor of 
fv « 0.95 (< 10“®-® Zq) and fv « 0.69 (< 10“® ® Zq). The 
volume filling fraction, fv, is defined as 


_ ^-^subset tttz /ftj 


(9) 


For gas that is sufficiently metal-enriched to form Pop II stars, the 
volume filling fraction keeps increasing and reaches fv ~ 0.05 
(> 10“®'® Zq) and fv « 0.29 (> 10“®'® Zq) at 2 « 10.5. 

The upper panel of Fig. 6 shows the volume filling fraction of 
the gas as a function of metallicity at four different characteristic 
times. Overall, the metal-enriched volume fraction increases over 
time, peaking at Z « 2 — 3 x 10 -"Zq and declining at high metal- 
licities Z > 10“® Zq. We also show the average metallicity as a 
function of gas density at the four selected times. The metallicity 
distribution sensitively depends on the gas density at the time when 
the SN is about to explode. The initial conditions for a SN explosion 
are set by the photoheating from its progenitor star. If the stellar ra¬ 
diative feedback is strong enough to evacuate the gas around the 
star, lowering the gas density to nn = 0.1 — 10“® cm“®, this low 
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Figure 7. Assembly history of the emerging target galaxy. Clockwise from upper left. Mass growth of the target galaxy, growth rate of various halo components, 
H I/H II fractions, and the average gas temperature. Starting from ~ 5 X 10® Mq at 2 : ^ 28, the galaxy grows in mass and reaches the threshold for 
atomic hydrogen cooling at 2 Ri 11.7. The final virial mass at the end of the simulation, 2 ^ 10.5, is Mvir ~ 10® Mq. The stellar content of the galaxy 
is dominated by Pop II, with a total mass of M* Ri 6 X 10® Mq, formed with an effective star formation rate of ~ 10“^ Mq yr“®. Through a complex 
inteiplay between metal ejection and re-accretion along the filaments of the cosmic web, roughly ~ 40 Mq in metals are present inside the galaxy at the end 
of the simulation. 


density gas is first likely to be polluted. For lower mass progenitor 
stars, on the other hand, due to their weaker ability of evacuating 
gas in their surroundings, such high density gas is preferentially 
enriched by the explosion. 


This trend is clearly displayed in the bottom panel of Fig. 6. 
The progenitor mass of the first SN, triggered at fn = 133 Myr, 
is 25 Mq. The outflow, driven by the photoheating of the 25 Mq 
Pop III star, is not strong enough to evacuate the surrounding dense 
gas. Therefore, the subsequent SN ejecta encounter this high den¬ 
sity gas first, giving rise to the metallicity peak at wh = 100 cm“® 
(see solid line in the bottom panel of Fig. 6). On the contrary, suffi¬ 
ciently strong radiative feedback from massive progenitor stars sets 
the stage for the following SN explosions to occur in low-density 
gas, with riH = 0.01 — 0.1 cm“®. This low density gas is then pref¬ 
erentially enriched, corresponding to the cases shown as dotted and 
dashed lines in Fig. 6. In summary, the derived trend is as follows: 
a SN with a low mass Pop III progenitor tends to pollute dense gas 
while a massive progenitor followed by a SN is more likely to en¬ 
rich low density gas. This trend, however, disappears as both low 
and high density regions are equally polluted as SN events become 
more frequent. 


3.2 Assembly of the first galaxy 

3.2.1 Mass evolution 

The first star in the Lagrangian volume of the emerging target 
galaxy forms at 2 « 28 out of the primordial gas within the main 
progenitor minihalo of mass ~ 5 x 10® Mq. After its short life¬ 
time of 2.7 Myr, the star collapses into a BH, acting as a HMXB 
source for the next 2 Myr. Stellar and BH feedback evacuate the 
gas by reducing the baryon fraction of the halo from /bar = 0.12 at 
2 « 28 to essentially zero at 2 « 25. Afterwards, the baryon frac¬ 
tion again increases, and at 2 « 16 gas densities have increased 
to tih = 10^ cm“®, thus enabling further star formation. Subse¬ 
quently, the baryon fraction remains above /bar = 0.05 throughout 
since 2 < 16, even in the presence of ongoing star formation, as 
the halo is sufficiently massive, Mvir ^ 5 x 10® Mq, to hold onto 
its gas despite feedback. 

The first Pop II stars in the target halo form at 2 = 13.1 out 
of gas that was enriched by a 33 Mq Pop III star that ended its 
live as a CCSN, thus initiating the prolonged chemical enrichment 
process. Until 2 = 13.3, the halo grows in mass mainly via smooth 
accretion along the filaments of the cosmic web and through minor 
mergers. At 2 ~ 13.3, the halo experiences a major merger and 
doubles its mass to 2 x 10^ Mq at 2 « 12.6, over a time span 
of ~ 25 Myr. The average growth rate during this period is M ~ 
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Figure 8. Morphology in the target galaxy at 2 : ps 10.5. Clockwise from upper left: Dark matter overdensity, hydrogen number density, gas metallicity, and gas 
temperature, averaged along the line of sight within the central ~ 80 kpc (comoving). Dashed white circles denote the virial radius, which is rvir ~ 1-1 kpc 
at this time, with a corresponding virial mass of Mvir ~ 8 x 10^ Mq. The positions of the Pop III remnants are marked as filled circles in the upper-right 
panel: BHs as yellow, neutron stars left behind by lower-mass CCSNe as blue, whereas the location of Pop II clusters are shown in white. Due to the ongoing 
star formation inside the galaxy, the central gas is substantially evacuated into the void that is roughly perpendicular to the filaments of the cosmic web, 
resulting in gas densities nn $ 10“^^ cm“^ with temperatures well above 10^ K due to SN blastwaves. The central region is substantially metal-enriched to 
and all gas inside the virial radius is polluted above Z = 10 ® Zq . 


O.5M0yr“^. At « = 11.7, the halo exceeds 3.5 x ItF Mq, and 
fulfills the criterion for the onset of atomic hydrogen cooling (see 
Fig. 7). At the end of the simulation at 2 « 10.5, the target halo has 
reached a virial mass of ~ 8 x 10^ Mq, after undergoing another 
major merger at a = 11.3. 

Fig. 7 also shows the growth rates onto the target halo, the 
evolution of the average gas temperature, and the FI I/H II frac¬ 
tion. The gas is photoheated by Pop III and Pop II stars to tempera¬ 
tures of ~ 10^ K, temporarily reaching even higher temperatures in 
the wake of SN explosions or due to X-ray heating from accreting 
BHs and HMXBs. Early on, the gas can quickly recombine once an 
emission source is turned off. However, at < 2 ~ 12, about 50% 
of the gas remains ionized by continuous star formation activity. 
The morphology of the emerging galaxy at 2 « 10.5 is shown in 
Fig. 8. Along a direction roughly perpendicular to the filaments of 
the cosmic web, the gas is highly evacuated to nn < 10“^^ cm“® 
by Pop II stellar feedback, heating it to T ~ 10® — 10^ K and 
enriching it to 2 ~ 10“^ ® Zq (see the bottom panels of Fig. 7). 


3.2.2 Star formation history 

Until 2 = 12.5, star formation in the target halo is dominated by 
Pop III stars, but afterwards Pop II stars form out of metal-enriched 
gas, and become the dominant stellar population. By the end of the 
simulation, the emerging target halo has hosted 15 accreting BHs, 
2 HMXBs, 40 Pop III CCSNe, possibly leaving neutron stars be¬ 
hind, and 31 Pop II clusters. The total mass of the stellar remnants 
increases from ~ 100 Mq at 2 « 28 to 3 x 10"^ Mq at 2 « 10.5 
(see the left-upper panel of Fig. 7). However, not all of them are 
formed inside the halo. In order to determine their birth places, we 
trace their trajectories in Fig. 9, by recording their distances from 
the center of the target halo. We find that 30% of the remnants, cor¬ 
responding to 15% in mass, formed in neighboring haloes and were 
subsequently accreted by the target halo. 

One of the key physical quantities of the galaxy is the total 
stellar mass, needed to predict its observational signature. Given 
the extremely short lifetime of Pop III stars, of order a few Myr, 
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Figure 9. Trajectories of all stellar remnants formed during the assembly 
process of the target galaxy. The dotted line denotes the virial radius of the 
halo. As can be seen, about 30% of the remnants, corresponding to 15% 
of total remnant mass, formed in neighboring haloes and were incorporated 
into the target halo. The majority of Pop II clusters are bom in-situ, out of 
gas inside the galaxy, and only 3 clusters arise from mergers. 



10“* 0.01 1 100 lO"* 


we only consider Pop II clusters to estimate the stellar mass, com¬ 
puted as M*,pop II oc rn °‘m dm, where the IMF slope is 

’ ^min 

a = 1.35, and the lower mass limit mmin = 0.1 Mq. The upper 
mass limit, on the other hand, is time variable and depends 

on the age of the Pop II cluster, such that only stars with a main- 
sequence lifetime larger than the age of the cluster contribute to 
the persistent stellar mass. The estimated stellar mass in the halo 
starts from 540 Mq at 2 « 13 and increases to 6 x 10^ Mq at 
2 ; « 10.5. The corresponding star formation rate (SFR) is a few 
~10-"MQyr-\ 

The small characteristic mass for Pop III stars, Mchar = 
20 Mq, adopted in this simulation, results in a qualitatively very 
different star formation and galaxy assembly history, compared 
with the Mchar = 100 Mq, which was often used in previous 
first galaxy formation work (Bromm & Yoshida 2011). For in¬ 
stance, our simulation employs the same initial conditions as those 
in Greif et al. (2010) and Jeon et al. (2012), both of which used 
Mchar = 100 Mq for Pop III stars and neglected the Pop III/Pop II 
transition. Greif et al. (2010) showed that the strong feedback from 
a single PISN and its progenitor, formed at the beginning of their 
simulation at 2 ~ 28, evacuates most of the gas out of the host 
halo. Except for one additional star at 2 « 14, the halo remained 
without further star formation for the next ~ 300 Myr after the 
PISN. More stars formed in the simulation of Jeon et al. (2012), 
where the target halo hosted a total of 8 stars in the presence of 
the positive X-ray feedback from an accreting BH source. To the 
contrary, the same target halo here forms a total of 52 Pop III stars 
and 31 Pop II clusters by the end of the simulation. This striking 
difference in the star formation history reflects the relative strength 
of radiative and mechanical feedback from Pop III stars, depending 
on their masses: photoheating from massive Pop III stars strongly 
evacuates the gas and suppresses further star formation, as demon¬ 
strated in Greif et al. (2010) and Jeon et al. (2012). 


Figure 10. Top'. Density-temperature phase diagram of the gas inside the 
target galaxy at 2 Ri 10.5. The color coding indicates gas metallicity. 
Thanks to the metal cooling, high density gas (nn > 100 cm“®) is able 
to cool down to the CMB temperature floor at ~ 30 K. Bottom: Density- 
metallicity phase diagram. The color coding now indicates gas temperature. 
The metallicity spread in the range Z = 10“^ — 10“^ Zq for high den¬ 
sity gas, which is mainly responsible for subsequent star formation, implies 
that both Pop III and Pop II might form in the enriched gas, but this de¬ 
pends on the choice of critical metallicity (see main text). Note that the gas 
with nil < 0.01 cm“® is highly enriched to ~ 0.01 Zq, due to the re¬ 
cent SN explosion in the rarified region established by photoheating from 
its progenitor. 


3.2.3 Metal enrichment 

The 40 Pop III and 248 Pop II CCSNe that exploded in the target 
halo produce a total of 28 Mq and 70 Mq in metals, respectively. 
Most (60%) of the metals are permanently ejected from the target 
system into the diffuse IGM, leaving only 40 Mq of metals in the 
galaxy, supplied partially by the re-accretion of metals, or provided 
by in-situ SN events. The mass-weighted, average gas metallicity 
in the galaxy is <Z>= 7 x lO”'^ Zq at 2 Ri 10.5, whereas that of 
the densest gas, out of which new stars form, is higher by a factor 
of~ 30. 

A notable difference from many previous studies is that all the 
heavy chemical elements inside the target galaxy are produced by 
CCSNe, and not in PISN events. Despite the smaller yield y = 0.05 
for CCSNe, in contrast io y = 0.5 for PISNe, and the smaller 
masses of Mchar = 20 Mq for Pop III stars, the total mass in met¬ 
als, now provided by a large number of Pop III CCSNe, is compara¬ 
ble to the metal production by a single 200 Mq PISN. The result is 
a level of overall metal enrichment similar to previous studies. For 
instance, the central gas within a Mvir ~ 10® Mq halo at 2 « 10 
is enriched to an average of ~ 10“® Zq by a single 200 Mq PISN 
(e.g. Greif et al. 2010). With a series of PISN events, the central gas 
cloud in Mvir = 5 x 10^ — 10® Mq haloes could be enriched to 
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Figure 11. Enriched gas mass fraction within the virial radius of the galaxy 
as a function of metallicity at 2 ^ 10.5. It can be seen that the peak remains 
roughly at Z = 10“'* — 10“^ Zq for all densities, except for the lowest 
density gas, with npi < 0.01 cm“^, which shows a peak at > 10“® Zq. 
This is due to a recent SN explosion, which is preferentially propagating 
along the underdense region of the galaxy. 

[Z/H] = —2.5 — 3.0, as well (e.g. Wise et al. 2012; Biffi & Maio 
2013). Therefore, we emphasize that while the choice of charac¬ 
teristic Pop III mass results in substantially different star formation 
histories, as discussed above, the degree of metal enrichment might 
be similar whether a first galaxy is enriched by a single or a few 
energetic PISNe, or by numerous CCSNe. 

To gain further insight into the stellar population mix in our 
target galaxy, in Fig. 10 we show the gas temperature (top) and 
gas metallicities (bottom) inside the galaxy as a function of hy¬ 
drogen density at 2 « 10.5. The metallicities of the high-density 
gas (riH > 10 ® cm“®), which is the main reservoir for star forma¬ 
tion, are in the range 10~^ Zq to 10“® Zq. As a result, instead of 
having a clear termination of Pop III star formation, both popula¬ 
tions may coexist within the galaxy. Lowering the critical metal¬ 
licity, however, as suggested by the dust-cooling theories, would 
lead to exclusively Pop II star formation inside the target system. 
Gas temperatures reach the floor set by the CMB at ~ 30 K (hori¬ 
zontal line in the top panel of Fig. 10), due to efficient metal fine- 
structure line cooling. Without it, the gas would only be able to cool 
to ~ 200 K via molecular hydrogen cooling. Overall, 40% of the 
gas mass within the galaxy is enriched above the critical metallicity 
for fine-structure cooling, and no primordial gas (Z < 10“® Z©) 
remains. As shown in Fig. 11, across all densities, the gas metal¬ 
licity peaks at 10“^ — 10“® Zq, excluding the lowest densities, 
n-H < 0.01 cm“®, where the peak is located above 10“® Zq. This 
difference is associated with recent metal pollution in low density 
regions, which are already highly evacuated by photoheating from 
the massive SN progenitor. 

The bright afterglow spectra of possible Pop III gamma-ray 
bursts (GRBs), with absorption lines imprinted by intervening 
metal-enriched clouds, would be an ideal tool for probing the metal 
enrichment of the early Universe simulated here (e.g. Wang et al. 
2012; Macpherson et al. 2013). Specifically, the equivalent widths 
of the metal absorption spectra may allow us to distinguish whether 
the first heavy elements were produced by a PISN or a CCSN, and 


in turn constrain the Pop III IMF (e.g. Wang et al. 2012). As we 
have discussed in Section 3.1, the diffuse, low column density, IGM 
may be the preferred hunting ground for a possible PISN enrich¬ 
ment signature, whereas the higher-column density systems may 
be dominated by the conventional CCSN signature. 

3.3 Stellar and galactic archeology 

A promising method to constrain the properties of primordial stars 
is to probe the metal-poor stars in the halo of our Milky Way, or in 
nearby dwarf satellite galaxies. This approach, often termed stel¬ 
lar or galactic archeology (e.g. Frebel et al. 2010; Brown et al. 
2012; Karlsson et al. 2013; Vargas et al. 2013; Milosavljevic & 
Bromm 2014), is based on the assumption that such metal-poor 
halo stars were bom out of gas that was polluted by a single or 
a few Pop III SNe, and have survived until the present day, thus 
preserving the signature of the first generation of stars (e.g. Simon 
et al. 2015). In particular, investigating the metallicity distribution 
of metal-poor stars, the so-called metallicity distribution function 
(MDF) (see Beers & Christlieb 2005 for a thorough review), places 
constraints on several key physical quantities, such as the critical 
metallicity, Zcrit, and the characteristic mass of Pop III stars. The 
MDF for the metal-poor Milky Way halo stars exhibits a peak at 
[Fe/H] « —2, and rapidly declines toward lower metallicities with 
a sharp cutoff at [Fe/H] « —4. 

The semi-analytical study by Salvador! et al. (2007), where 
a Monte Carlo methodology based on an analytic Press-Schechter 
merger-tree is employed, shows that their model for which Zcrit = 
10“"^ Zq and Mpop in,char = 200 Mq can nicely reproduce the 
observed Galactic halo MDF, except for the existence of hyper 
metal-poor stars below [Fe/H] = —5. Fig. 12 shows the MDF 
of Pop II stars in our simulation, within the target halo (red) and 
within the more extended high-resolution region (black), evaluated 
at the end of our simulation at 2 « 10.5. We should mention 
that stellar metallicities, Zs, are different from the gas metallicities 
discussed in previous sections. The stellar metallicity is inherited 
from the gas metallicity when a star is bom out of its birth cloud. 
For the conversion between Z^ and [Fe/H], we use the relation 
[Fe/H] ~ log,o(Zc/Z 0 )-O .25, adopting the relative iron frac¬ 
tion produced by Pop III stars (Heger & Woosley 2010). The cut¬ 
off at [Fe/H] = —3.75 arises from the critical metallicity Zcrit = 
10 -®'® Zq adopted here, below which Pop II stars cannot form. 
Thus, all stellar metallicities exceed [Fe/H] = —3.75 by design, 
and the maximum is at [Fe/H] = —3.75 and [Fe/H] = —3.55 for 
the target halo and the high-resolution region, respectively. 

The stellar metallicity is determined by a complex interplay 
between the strength of photoheating by Pop III stars and the ability 
of haloes to retain their gas. In a halo with a shallow potential or in 
a halo hosting a massive Pop III star, the stellar metallicity of Pop II 
stars tends toward the metallicity of the reincorporated gas that was 
ejected out of the host halo by a previous SN blastwave. On the 
other hand, if photoheating is not strong, as is usually the case for 
the low mass Pop III stars (Mpop iii < 40 Mq), or if a host halo 
is massive. Pop II stars are likely to be bom out of the gas that is 
highly enriched, because after a SN explosion the ejected metals 
remain in the central region and the central gas density remains 
high, promptly triggering Pop II star formation. 

We should caution that our MDF is not directly comparable 
to the observed MDF, particularly to the low-metallicity tail of the 
observed MDF ([Fe/H] < —3.75). The main reason for the lack of 
these extremely metal-poor (EMP) stars in the simulated MDF is 
due to the critical metallicity we adopt. Hence, in order to include 
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such stars we might use a lower value of the critical metallicity. 
The other way of explaining EMP stars would be to assume that 
low mass, long-lived metal-free stars (Mpop m < 0.8 Mq) were 
formed, and that their surfaces were polluted externally later on by 
accretion of gas enriched with metals (e.g. Shigeyama et al. 2003; 
Trend & Shull 2010; Johnson & Khochfar 2011). 

One of the interesting open questions in high-z galaxies is 
whether or not we can find any fossils of the first galaxies in the 
Local Group, especially in local dwarf galaxies (see Tolstoy et al. 
2009 for a review; see also e.g. Bovill & Ricotti 2009; Simpson 
et al. 2013; Milosavljevic & Bromm 2014). One widely discussed 
scenario for the origin of dwarf galaxies invokes the suppression 
of their gas supply due to reionization (e.g. Salvadori & Ferrara 
2009). Specifically, most stars are argued to have formed before 
reionization in haloes with Mvir = 10® — 10® Mq, giving rise to 
systems that have then passively evolved into a subset of the dwarf 
galaxies. Next to reionization, photoheating and SN feedback act 
to evacuate the gas from the haloes. This scenario would naturally 
explain both the old age of the extant stellar populations, and the 
absence of any ongoing star formation in dwarf galaxies. However, 
the probability of a galaxy formed near the reionization era to re¬ 
main intact and isolated, not undergoing mergers, is small, of order 
~ 10% (Sasaki 1994). Therefore, only a subset of the present-day 
Local Group dwarfs might be suitable candidates for comparison 
with the first galaxies (Frebel & Bromm 2012; Frebel et al. 2014). 

Fig. 13 shows the relation between the mean metallicity of 
metal-poor stars and their masses in dwarf galaxies of the Milky 
Way and M 31 (Kirby et al. 2013). Assuming that stellar metallic- 
ities remain unchanged until today, the estimated mean metallic¬ 
ity of our simulated galaxy is <[Fe/H]>~ —3.3, which is lower 
by « 0.5 dex than what is expected from the empirical rela¬ 
tion, <[Fe/H]>Ri —1.69+ 0.30 log^Q(Mt/10® Mq) (see equ. 4 in 
Kirby et al. 2013). However, our simulation ends before reioniza¬ 
tion is complete. Indeed, further star formation and metal enrich¬ 
ment will likely occur, leading to Pop II stars with higher metal- 
licities. The effect of this missing activity may be to raise the av¬ 
erage metallicity, possibly up to levels consistent with the empiri¬ 
cal Kirby-relation. Therefore, tracing the evolution of the galaxy to 
2 ; ~ 5, and ideally even beyond that, might be necessary to directly 
compare simulation results with observations. 

To infer the star formation history of dwarf galaxies, Web¬ 
ster et al. (2014) have recently performed simulations that fol¬ 
low the evolution of an isolated system with a dark matter mass 
of ~ 10^ Mq for 600 Myr, accounting for metal enrichment by 
Type II and Type la SNe. They found that their simulated system 
has a mean metallicity of [Fe/H] = —2.1 and a stellar mass of 
2 X 10® Mq after 600 Myr, suggesting Segue 1 as the closest ob¬ 
served match (see also Frebel et al. 2014). However, their simula¬ 
tion is performed in a small region (< 400 pc) with a rather idealized 
setup. In particular, their high mean metallicity is caused by the as¬ 
sumption of an enrichment floor that is linearly increasing with the 
number of SNe. On the other hand. Pop II SN enrichment in our 
simulation does not impose any artificial metallicity floor, thus en¬ 
abling lower metallicity conditions. 

3.4 Luminosity 

One of the imminent challenges for the next generation of tele¬ 
scopes, such as the JWST, is to detect the first light from high- 
redshift dwarf galaxies in the early Universe. In order to assess the 
detectability of these sources, we compute the expected bolomet- 
ric and line emission luminosities for the simulated galaxy. Due to 
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Figure 12. Metallicity distribution function (MDF) of Pop II stai's formed 
by the end of the simulation. The MDF peaks at Zs = Zq and 

Zs = 10 “^'^ Z 0 , respectively, for the stars inside the target halo (red) and 
for all Pop II stars in the high-resolution region {black). The lower cutoff 
is a consequence of the adopted critical metallicity, which may be lower in 
reality. 



Figure 13. Stellar mass-mean metallicity relation from observations 
of dwarf galaxies around our Milky Way and M 31. They are 
strongly correlated, well fitted by the formula <[Fe/H]> 5 ^ —1.69 + 
0.30logjQ(Af*/10® Mq), shown as the dashed line (Kirby et al. 2013). 
For comparison, we also plot the estimate for the mean metallicity of our 
simulated galaxy {red dot). This value is lower than the empirical relation 
by 5 ^ 1 dex, implying that the simulated system cannot give rise to the 
observed local dwarfs. However, see the main text for further discussion. 

the short lifetime of Pop III stars, we only consider emission from 
Pop II stars. Assuming an instantaneous starburst, the bolometric 
luminosity from a Pop II cluster is represented by the following 
equation, 

/’TTlmax 

Lfioi.Pop II OC / m~°‘(j>{m) dm, (10) 

where (f>{m) is the mass-luminosity relation. Depending on the 
mass of a star, a different power-law is used: L oc m®'® (m < 
0.43 Mq), L oc m* (0.43 Mq < m < 2Mq), L oc m® ® 
(2 Mq < m < 20 Mq), and L oc m (m > 20 Mq). Additionally, 
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Figure 14. Bolometric luminosity {left) and recombination line luminosity from Ha and Lya (right) for haloes hosting Pop II clusters. The target galaxy 
has a bolometric luminosity of LboI = 3.5 x 10^ Lq. Except for this tai'get galaxy, the bolometric luminosities of the other haloes are of order a few times 
lO^Lo. The numbers next to the red dots in the left panel indicate the mass-to-light ratios, which are quite large, similar to that of the local dwarf galaxies. 
The recombination line luminosities from Ho and Lyo are smaller than LboI by 1 — 4 orders of magnitude. 


we compute the nebular emission, i.e. the luminosity in recombina¬ 
tion lines, particularly for Ha and Lya, from the surrounding gas 
that is ionized by star formation (e.g. Johnson et al. 2009), 
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where j is the temperature-dependent emission coefficient of the 
lines (Osterbrock & Ferland 2006), and /e and /hii are the frac¬ 
tions of electrons and HII ions of the SPH particles inside the virial 
volume of a given halo. 

Fig. 14 shows that the most massive halo with a mass Mvir ~ 
8 X 10^ Mq, hosting the target galaxy, has a bolometric luminosity 
of Z/Boi = 3.5 X10® Lq . We note that 20% of Pop II clusters (out of 
34) are younger than 20 Myr, and these are the main contributors to 
the luminosity. For the other less massive haloes in the mass range 
of Mvir = [1-2 X 10® Mq, 3 X 10^ Mq], the bolometric luminosity 
varies in the range of LboI = [8.7 x 10^ Lq, 2.6 x 10® Lq]. Our 
estimates are consistent with those in Wise et al. (2014), who use 
a more sophisticated stellar population model to derive bolometric 
luminosity. As stated above, we only consider the luminosity from 
Pop II clusters. We should note, however, that the simulated target 
galaxy continues to form Pop III stars as well, out of the gas with 
metallicity below Zcrit = 10“®'® Zq. At z ~ 10.5, the galaxy 
hosts 6 Pop III stars with a total luminosity of 1.4 x 10® Lq, corre¬ 
sponding to 40% of the Pop II bolometric luminosity. Thus, adding 
the contribution from Pop III stars, the combined luminosity may 
increase roughly by a factor of ~ 1.4. 

Using our bolometric luminosity estimates, we can predict the 
observed flux at a rest-frame wavelength of Ae = 1500A, char¬ 
acteristic of the soft (non-ionizing) UV continuum. The observed 
continuum flux for an unresolved object is given by (e.g. Oh 1999; 
Bromm et al. 2001) 


/4Ao) 
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where dr is the luminosity distance to a source at redshift z, which 
is ~ 100 Gpc at 2 = 10. The predicted flux for the galaxy 
(Mvir = 8 X 10^ Mq) in our simulation is ~ 10“® njy for 
LboI = 3.5 X 10® Lq. This is smaller by two orders of magnitude 
than the detection limit of the Near Infrared Camera (NIRCam) on¬ 
board the JWST for extremely deep imaging (e.g. Gardner 2006). 
Therefore, more promising targets for detection with the JWST will 
be more massive, luminous, systems with Mvir = 10® — 10®® Mq 
(e.g. Pawlik & Schaye 2011). The other possible way of detecting 
those small galaxies is through gravitational lensing with a boost in 
brightness by a factor of /r ~ 100 (e.g. Zackrisson et al. 2012, 
2015). However, magnifications of ^ > 100 are extrmely rare, 
mainly due to the low probability of having sight lines for the grav¬ 
itational lensing, P(> n) oc (e.g. Pei 1993), requiring very 
large survey areas to find the lensed galaxies. 

On the other hand, the nebular emission in recombination lines 
are Lna = 4 x 10® Lq and I/Lya = 4.8 x 10® Lq, respectively, 
both of which are even weaker than the bolometric luminosity by 
2 — 3 orders of magnitude. This is mainly due to the low gas den¬ 
sities as a result of the photoheating from the newly forming stars, 
resulting in low recombination rates. In our estimate, we have as¬ 
sumed escape fractions for Ha and Lya photons of unity, rendering 
the calculated line luminosities upper limits. Although the recom¬ 
bination luminosity in Lya is higher than that in Ha, in reality 
Lya photons are more likely to scatter in the neutral IGM (e.g. Di- 
jkstra & Wyithe 2007; Smith et al. 2015), and thus the observed 
flux in Lya may be substantially lower than implied by the in¬ 
trinsic luminosity. Because of the hard spectra of massive Pop III 
stars, strong nebular emission, particularly in the He II A1640 line, 
has been suggested as a possible observable signature of metal-free 
Pop III galaxies (e.g. Schaerer 2003; Dawson et al. 2004; Nagao 
et al. 2008; Johnson et al. 2009; Cai et al. 2015). However, as im¬ 
plied by our simulation, the gas inside the first galaxies is likely 
already metal-enriched, and thus their stellar content is expected to 
be dominated by Pop II, instead of Pop III, stars, making it more 
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difficult to use the recombination lines as a probe to detect such 
metal-enriched small systems. 

Finally, a complementary method is to look for counterparts of 
the first galaxies in the Local Group. Assuming, probably unrealis¬ 
tically so, that star formation in our simulated galaxy is completely 
quenched after 2 : « 10.5, and that it has passively evolved until the 
present day, the total stellar mass is 1.4 x 10® Mq at 2 = 0 for the 
IMF we use. This is comparable to the stellar mass M* = 900 Mq 
of Segue 2, the least luminous dwarf galaxy known in the local 
Universe (Kirby et al. 2013). More high-resolution observations of 
the faintest galaxies are needed to further elucidate the connection 
between the fossils of the first galaxies and local dwarf galaxies. 


4 SUMMARY AND DISCUSSION 

We have carried out a radiation hydrodynamical simulation of a 
galaxy reaching a mass of Mvir = 8 x 10^ Mq at 2 « 10. This 
simulation tracked ionizing photons from individual Pop III stars 
and from Pop II clusters. High spatial resolution allowed us to 
explicitly account for the mechanical and chemical feedback from 
CCSNe and PISNe. Additionally, we followed the transport of 
X-ray photons from accreting black holes (BHs) and high mass 
X-ray binaries (HMXBs). Motivated by recent works on primordial 
star formation, which suggest Pop III stars with a few tens of 
Mq, we adopted a characteristic mass, Mchar = 20 Mq for the 
Pop III IMF. This is substantially smaller than Mchar = 100 Mq, 
which has been used in most previous studies. We also investigated 
the metal enrichment history in the early Universe, including the 
transition in star formation mode from metal-free Pop III stars to 
low-mass, metal-enriched Pop II stars at the critical metallicity 
.^crit = 10-®'® Zq. Our main results are as follows. 

(i) Adopting a moderate characteristic mass Mchar = 20 Mq, 
instead of Mchar = 100 Mq, for the Pop III IMF, results in a 
substantially different star formation history during the assembly 
process of the target halo that hosts the first galaxy. This is caused 
by the relatively weak stellar feedback, particularly photoheating, 
from the less massive Pop III stars. 

(ii) The gas inside the galaxy is enriched to Zavg ~ 10“® Zq 
by a total of 40 Pop III and 248 Pop II core-collapse supernovae 
(CCSNs). Despite the striking difference in the star formation 
history, the degree of metal enrichment by CCSNe is comparable 
to that by a single 200 Mq Pop III pair-instability supernova 
(PISN). 

(iii) The metals produced by two PISNe outside the main 
galaxy are mainly responsible for the metal enrichment in the dif¬ 
fuse IGM, setting a metallicity floor of Z = 10“"^ Zq. Therefore, 
signatures of PISNe in the early Universe may be found by studying 
the state of the IGM, e.g. employing bright Pop III CRB afterglows. 

(iv) There is no gas with metallicity Z < 10“® Zq inside the 
main galaxy. Given the critical metallicity Xcrit = 10“®'® Mq 
adopted here. Pop III stars continue to form out of the gas in the 
metallicity range 10“® Zq < Z < 10“®'® Zq. Thus, both Pop III 
and Pop II stars coexist. This conclusion, however, is uncertain, as 
the critical metallicity may well be lower, which would not allow 
any extended Pop III star formation. 

(v) Due to the rarity of HMXBs and their short life¬ 


times, in combination with the low BH accretion rates 
Mbh = 10“® — 10“®^MQyr“®, the impact of X-rays on 
the gas inside haloes that host X-ray sources and on the IGM is 
negligible in this simulation. 

(vi) The bolometric luminosity of the simulated galaxy is 
LboI ~ 3.5 X 10® Lq and the total stellar mass is 6 x 10® Mq 
at 2 « 10.5. Therefore, the UV flux of the galaxy, redshifted into 
the near-IR is 1.4 x 10“® njy, significantly less than the detection 
limit for NIRCam onboard the JWST. 

(vii) The estimated metal distribution function (MDF) of 
Pop II stars in the galaxy peaks at [Fe/H] = —3.75 and the 
stellar metallicities are in the range [Fe/H] = [—3.75, —2.7] 
with a mean metallicity of [Fe/H] = —3.3. The stellar mass of 
M* = 1.4 X 10® Mq of the simulated galaxy, passively evolved 
to the present day, is comparable to that of Segue 2, the least 
luminous dwarf known in the local Universe. 

It has often been pointed out that including a physically moti¬ 
vated description of stellar feedback is crucial for modeling galaxy 
formation, especially for small galaxies, which are susceptible to 
stellar feedback due to their shallow potential wells (e.g. Wise et al. 
2012; Kim et al. 2013; Hopkins et al. 2014). We confirm that stel¬ 
lar feedback plays a crucial role in shaping the properties of the 
first galaxies: by reducing the baryon fraction (/bar ~ 0.05), cor¬ 
responding to 15% of the cosmic mean, and by enriching the gas to 
an average of <Z> ~ 10“® Zq, making Pop II star formation the 
dominant mode inside the first galaxies. However, galaxies may 
continue to form Pop III stars because not all gas is enriched above 
the critical metallicity. 

In our simulation, we needed to fix a large number of numer¬ 
ical input parameters, such as the random sampling of the Pop III 
stellar masses from the underlying IMF. Given that computational 
resources limit us to only one such realization from a very large 
ensemble of likely parameter choices, the question arises what as¬ 
pects of our simulation are contingent on the specific choice, and 
which results are likely robust. An example of a robust feature is 
that the first galaxies were already metal-enriched, such that their 
stellar content is expected to be Pop II; the precise level of this 
‘bedrock metallicity’, on the other hand, is uncertain, but likely 
of order Zavg > 10 “'*Zq at 2 « 10. An example of a contin¬ 
gent, numerical-parameter dependent feature is the impact of X-ray 
feedback from early HMXB sources. We mention that X-ray feed¬ 
back is negligible in our simulation, but the population of HMXBs, 
which are likely the strongest X-ray sources, and their luminosity 
and duty-cycle are unknown at high redshifts. We refer the reader 
to Jeon et al. (2012) for discussion of a case where X-ray feedback 
is more significant. 

We show that the stellar mass fraction in halos is /* = 
M*/Mvir ~ 10“^, which is smaller by an oder of magnitude 
than the values found in other studies, for example Behroozi et al. 
(2013) and O’Shea et al. (2015), using abundance-matching tech¬ 
niques and grid-based simulations, respectively. However, the scat¬ 
ter from their work is large, such that /, ranges from 10“® to 10“® 
at a given mass of Mvir = 10® Mq. On the other hand, our esti¬ 
mated for ft is in a good agreement with the results from Garrison- 
Kimmel et al. (2014), in which they use an updated abundance¬ 
matching technique and reproduce the stellar mass functions of the 
Milky Way and Andromeda satellites more accurately. The fact that 
near the end of our simulation the baryon fraction of the target halo 
is significantly reduced to 15% of the cosmic mean, leading to a 
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reduction in star formation, robustly implies a highly feedback- 
regulated regime. Such inefficiency of star formation could be at¬ 
tributed to how we implement stellar feedback, especially the ef¬ 
ficiency of SN feedback. However, even if the stellar mass were 
increased by a factor of 10 by adopting a different parameter imple¬ 
mentation, this system would still be observationally out of reach 
for the JWST, at least in the absence of gravitational lensing. 

Whether fine-structure line cooling or dust-continuum cool¬ 
ing was mainly responsible for Pop II star formation is still de¬ 
bated (e.g. Omukai 2000; Bromm et al. 2001; Schneider et al. 2002; 
Bromm c& Loeb 2003b; Omukai et al. 2005; Schneider & Omukai 
2010; Schneider et al. 2012; Chiaki et al. 2014; Ji et al. 2014). If 
we had adopted a critical metallicity Zcrit,dust = 10“®'® Z©, set by 
dust continuum cooling, the assembly history of the target galaxy 
would have proceeded quite differently. Adopting a lower critical 
metallicity than Zcrit = 10“®'® Z© would boost the relative im¬ 
portance of Pop II star formation, and long-lived Pop II stars could 
hence contribute to the total stellar mass. The estimated MDF of 
Pop II stars would also be changed depending on the critical metal¬ 
licity. The nature and observational properties of the first galax¬ 
ies thus sensitively depend not only on the mass scale of the first 
stars, as has previously been pointed out (Bromm & Yoshida 2011), 
but also on additional astrophysical parameters, such as the critical 
metallicity. 

Considerable uncertainties remain in the treatment of metal 
mixing. We used a subgrid model (Greif et al. 2009) in which 
we assume that the efficiency of metal mixing below the resolved 
scales is determined by physical quantities at resolved scales. Shen 
et al. (2010) propose a new implementation to more accurately cap¬ 
ture the turbulent mixing in SPH. Further theoretical and obser¬ 
vational studies are required to better understand the metal mix¬ 
ing process. As discussed, exploring a large range of parameters 
is needed to derive representative characteristics of high -2 dwarf 
galaxies. However, performing such a comprehensive suit of high- 
resolution, radiation-hydrodynamical simulations is computation¬ 
ally challenging. 

As mentioned above, the direct detection of the simulated 
galaxy, with Mvir < 10® M©, will be out of reach even for 
the JWST, unless extreme, and probably unattainable, gravitational 
lensing magnifications are involved. Therefore, a more promising 
strategy to infer their nature may be to use “stellar" or “galactic" 
archeology by searching for metal-poor stars in the Milky Way or 
in nearby dwarf galaxies. To fully unleash its potential, both simu¬ 
lations and observational surveys need to be pushed further, a pro¬ 
cess that clearly has been set in motion now. 
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